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RESUME 


Nous  presentons  ici  une  nouvelle  methode  pour  obtenir  la  fonction 
de  coherence  mutuelle  d'une  impulsion  laser  se  propageant  au  milieu 
d 'aerosols.  A  partir  d’une  Equation  diffSrentielle  parabolique,  la 
fonction  de  coherence  mutuelle  est  dSterminSe  5  l'aide  d'une  transforma¬ 
tion  Fourier  'split-step'.  Cette  approche  nous  permet  d'Sviter  l'utili- 
sation  d' approximations  non  contr8l6es  comme  ce  fut  souvent  le  cas  dans 
les  etudes  prScedentes.  Des  algorithmes  numSriques  sous  forme  de  rela¬ 
tions  algSbriques  rScursives  permettent  d'Studier  le  ph€nom§ne  de  l'Slar- 
gissement  de  1' impulsion  caus£  par  la  diffusion  multiple.  (NC) 


ABSTRACT 

X 

A  new  method  of  obtaining  the  mutual  coherence  function  of  a  laser 
pulse  propagating  in  a  particulate  medium  is  described.  The  treatment, 
based  on  the  technique  of  the  split-step  Fourier  transform,  enables  us 
to  avoid  the  uncontrolled  approximations  often  used  to  solve  the  parabo¬ 
lic  differential  equation  of  the  mutual  coherence  function.  Simple 
numerical  algorithms  in  the  form  of  algebraic  recurrence  relations  are 
derived  for  the  study  of  the  phenomenon  of  pulse  broadening  induced  by 
multiple  scattering.  (U) 
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1.0  INTRODUCTION 

In  the  application  of  electro-optical  systems  such  as  laser 
rangers  and  laser  designators,  it  is  important  to  understand  the  propa¬ 
gation  characteristics  of  the  electromagnetic  pulses  emitted  from  the 
systems.  This  is  particularly  true  when  the  pulses  have  to  traverse 
a  turbid  medium  such  as  the  atmosphere  under  low-visibility  conditions 
produced  by  fog,  cloud  or  rain.  Apart  from  energy  attenuation,  the 
pulses  also  suffer  distortion  which  can  make  the  shape  of  a  received 
pulse  significantly  different  from  that  emitted.  In  addition,  in 
atmospheric  sensing  using  a  pulsed  optical  system  such  as  a  LIDAR 
device,  the  physical  information  of  the  medium  has  to  be  extracted  from 
the  analysis  of  the  characteristics  of  the  returned  pulses. 

In  the  study  of  the  pulse  propagation  characteristics  in  a 
medium  of  randomly  distributed  scatterers  one  may  consider  the  evolu¬ 
tion  of  the  mutual  coherence  function  as  in  the  work  of  Beran  (Refs.  1 , 
2)  and  Ishimaru  et  al  (Refs.  3-5) .  Alternatively,  the  technique  of 
temporal  moments  as  considered  by  Liu  and  Yeh  (Refs.  6-8) ,  or  the  direct 
solution  of  the  time  dependent  radiative  transfer  equation  (Refs.  9,  10) 
may  also  be  applied. 

In  the  present  study,  we  describe  a  new  method  of  deriving  the 
mutual  coherence  function  using  the  split-step  Fourier  method  (Refs.  11, 
12) .  Under  the  Markov  approximation  the  mutual  coherence  function 
satisfies  a  parabolic  differential  equation  (Ref.  13).  This  equation 
formally  resembles  a  time  dependent  Schroedinger  equation  with  an 
imaginary  'potential'  iP(p)  where  P(p)  is  the  Fourier  transform  of  the 
scattering  phase  function.  Under  various  different  approximations, 
this  equation  has  been  solved  analytically  by  Sreenivasiah,  Ishimaru, 
and  Hong  (Ref.  5)  as  well  as  by  Zardecki  and  Tam  (Ref.  14).  In 
parallel  with  these  efforts,  numerical  integrations  of  the  parabolic 
equations  have  also  been  attempted  (Refs.  5,  14). 
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Exploiting  a  technique  similar  to  the  time  evolution  operator 
expansion  used  in  solving  the  Schroedinger  equation  we  derive  a  recur¬ 
rence  relation  for  the  mutual  coherence  function  with  the  aid  of  the 
split-step  Fourier  method.  This  treatment  enables  us  to  avoid  both  the 
quadratic  expansion  of  P(p)  (Ref.  5)  or  the  Glauber  approximation  in 
Ref.  14.  In  other  words,  the  parabolic  differential  equation  for  the 
mutual  coherence  function  can  now  be  solved  without  invoking  simplifying 
assumptions . 

In  Section  2.0,  the  parabolic  differential  equation  of  the  mutual 
coherence  function  is  briefly  described.  The  new  method  of  solving  this 
equation  by  means  of  the  evolution  operator  expansion  and  the  split-step 
Fourier  scheme  is  discussed  in  Section  3.0.  In  Section  4.0,  the  problem 
of  stepwise  generation  of  the  mutual  coherence  function  is  considered. 
Numerical  results  for  laser  pulses  propagating  in  fog  are  presented 
in  Section  5.0.  The  work  was  performed  at  DREV  between  December  1979 
and  July  1980,  under  PCN  33A11,  Aerosol  Studies. 


2.0  PARABOLIC  EQUATION  OF  THE  MUTUAL  COHERENCE  FUNCTION 

A  convenient  mathematical  quantity  for  the  description  of  the 
electromagnetic  radiation  emitted  from  an  extended  source  is  the  mutual 
coherence  function  (Refs.l,  15).  In  a  random  turbulent  medium,  the 
propagation  of  the  mutual  coherence  function  has  been  shown  to  obey  a 
parabolic  equation  by  Beran  (Refs.l,  2).  The  validity  of  Beran's 
results  was  confirmed  by  Brown  (Ref.  16)  solving  a  Bethe-Salpeter  equa¬ 
tion  of  the  mutual  coherence  function.  For  a  medium  where  the  scatter¬ 
ing  is  due  to  randomly  distributed  particulate  matter,  a  similar  para¬ 
bolic  equation  applies  provided  the  scattering  events  are  predominantly 
in  the  forward  direction  (Refs.  3-5) . 
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Let  v(r,t)exp(-ia>ot)  be  the  output  field  at  a  receiver  located 
at  r  produced  by  an  input  field  u(t)exp(-iu)Qt)  at  the  transmitter 
located  at  the  origin.  The  input  and  output  fields  are  related  by  the 
equation  (Ref.  17) 


v(r,t)  = 


dui  U(ai)  H(o)  +  wo,  r,t)  exp(-iut) 


where  H(o),  r,t)  is  the  frequency  response  function  of  the  medium  and 
U(u>)  is  the  Fourier  transform  of  u(t).  The  two-frequency  mutual  coher¬ 
ence  function  r(u)j,a>2.*j,r2,t)  is  defined  as  the  ensemble  average: 


r(«1,w2»r1»r2»  O  =  <H *r^ »t)H* (w2 »^*2 *^) > ■ 


The  usefulness  of  the  mutual  coherence  function  lies  in  the  fact  that 
given  an  input  signal  whose  intensity  is  I^(t)  at  the  transmitter,  the 
output  signal  intensity  at  r  can  be  expressed  as  a  convolution 


01 

’  / 


dt’G(r,t-t*)I.(f) 


where  the  Green's  function  G(r,t)  is  related  to  r(Wj,u>2,rj,r2,t)  by  the 
Fourier  transform 


00 

G(r,t)  =  —  /  du,  r(aj,  +  u  ,  +  u>  ,r ,r ,t)exp(-iio.t) 

2*  J  d  !  o  2  o  d 


with  »  ujj  - 
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For  a  plane  wave  propagating  in  the  z  direction,  the  mutual 
coherence  function  at  two  points  r^  =  (p^  ,z)  and  on  the 

same  z  plane  will  depend  on  p d  =  and  2  but  will  be  independent 

of  tc  ®  j  (Pj  +  p^ •  Assuming  the  medium  to  be  time  invariant  and  the 
scatterings  to  be  uncorrelated  and  predominantly  in  the  forward  direction, 
it  can  be  shown  that  (Ref.  17)  the  mutual  coherence  function  written 
now  as  r(o)^,pj,z)  satisfies  the  following  parabolic  equation 


3_ 

3z 


ik 


2k 


v5- 


ikd  + 


at  -  p(0d) 


r(a)d,pd,z)  =  o 


[5] 


with  the  boundary  condition 


r(wd,Pd,z  =  0)  =  1 


[6] 


2 

In  eq.  5,  kQ  =  <o^/c ,  kd  =  u>d/c,  Vd  is  the  two-dimensional  Laplacian 
with  respect  to  pd,  o^.  is  the  volume  extinction  coefficient  and  P(Pd) 
is  the  Fourier  transform  of  the  scattering  phase  function  p(s): 


P(Pd)  - 


/ 


d  s  p(s) 


exp(-ikQS' 


pd} 


[7] 


with  os  denoting  the  volume  scattering  coefficient.  The  normalization 
condition 
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is  assumed. 

When  the  function  P(p)  is  quadratic  in  p,  an  analytic  solution  to 
eq.  5  is  available  (Ref.  5).  In  the  more  general  case,  an  analytical 
solution  has  also  been  obtained  using  Glauber's  eikonal  approximation 
(Ref.  -14).  Unfortunately,  Glauber's  approximation  does  not  reproduce 
the  time  asymmetry  in  the  pulse  shape  introduced  by  the  scattering. 


3.0  SPLIT-STEP  FOURIER  SOLUTION  OF  THE  PARABOLIC  EQUATION 
Let  us  introduce  a  mutual  coherence  function  y(a>,p,z)  defined  by 


y(u>,p,z)  -  r(w,p,z)  exp(atz  -  ikz) 


[9] 


where  we  have  simplified  the  notation  so  that  w  =  p  =  p^.  Equation 
5  becomes 


=  -iaA  +  P(p)J 


—  y(w,p,z)  =  |-iaA  +  P(p)  j  y(u,P,z) 
3z 


[10] 


2  2 

with  A  h  7^,  a  =  k/2kQ  and  the  boundary  condition 


Y(u),p,z  =  0)  =  1 


[11] 


<•'-*6 
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From  eqs.  3-4  the  function  y(w,p  =  0,z)  is  required  in  order  to 
evaluate  the  output  pulse  I(r,t). 

Consider  the  expansion 


y(u»,p,z  ♦  6z)  ■  y(“*Piz)  ♦  6z  +  —  — y  («*)  *  •••  1^2] 

9z  2!  3z 


By  means  of  eq.  10  one  can  write 

Y(w,ij,z+Sz)  ■  exp  [6z(-iaA  ♦  P(p))]  y(“»P>z)  t13l 


For  the  operator  exp  [<5z(-iaA  +  P(p))]  we  introduce  the  following 
approximation : 


exp[6z(-iaA  +  P(p))]  =  exp 


exp[-ia6zA]  exp 


[14] 


The  error  in  equating  the  two  sides  of  eq.  14  is  of  the  order  (6z)  . 

Let  g(p)  =  g(x,y)  be  an  arbitrary  square  integrable  function 
whose  Fourier  transform  is 


-00 


d2t  exp(-itp)  g(t)  =  F_1[g(t)] 


[IS] 
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Since 


oo 

Ang(p)  =  — ■ — r  ((  d2t  (-T2)n  exp(-it‘p)  g(t) 

(2*r 


[16] 


we  have 


exp(-iaSzA)  g(p)  =  — II exp  [  ia<S  zt2  ]  exp  ( -  i  t-p)g(t) 

C2*)2  r y 

=  F'1  [exp(ia6zA  )  F[g(p)]] 


[17] 


Combining  eqs.  13,  14  and  17  we  obtain  the  following  equation 


M  \ 

r 

-1 

y(w,pz,z+6z)  =  exp  [—  P(p)j 

F-1 

exp(ia6zt^)  F 

exp  (—  P(p)Jy(w,P,z) 

V  2  / 

- 

\  2  / 

[18] 


which  is  the  recurrence  relation  of  y(u>,p,z)  obtained  by  applying  the 
split-step  Fourier  algorithm  (Refs.  11,  12)  to  eq.  10. 


In  the  case  where  the  function  P(p)  possessed  cylindrical 
symmetry  about  the  z  axis,  the  following  simplification  can  be  made. 


y(u>,p,z)  JQ(Xp) 


[19] 
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where  J  (x)  is  the  Bessel  function  of  order  zero.  Equation  18  now 

O 

becomes 


y(<o,p,z+6z)  =  exp 


00 

J  dp  v  exp  ^p(p')^ 


2 

dX  exp(iaX  Sz)  XJQ(Xp) 


JQ(Xp  ')yO*>»p'  ,z) 


=  exp 


00 

^-P(p)^  dp,p’  exp^P(p')j  Y(«,P'.z) 


00 

/ 


dXX  exp(iaX25z)  JQ(Xp)  Jo(Xp') 


With  the  aid  of  the  formulas  (Ref.  18) 


/“  [(=) 


(aX26z) 


JQ(Xp)  J0(*P') 


2a6z 


/sin\  /p2  *  P  ' 2  \  j  /pp  1  \ 
\COS/  \  4afiz  /  °  y2aSz ) 


Eq.  21  can  be  reduced  to 
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In  the  case  of  a  laser  beam  passing  through  fog,  cloud  or  rain 
the  single  scattering  law  can  be  modelled  by  a  Gaussian  scattering 
phase  function  (Refs.  19,  20) 


P(s) 


°sa 


TF 


[23] 


where  t  is  the  projection  of  s  on  the  x-y  plane.  Corresponding  to  this, 
the  function  P(p)  defined  by  eq.  7  can  be  written  as 


P(p)  =  exp  (-k2  p2/4c0 


[24] 


To  simplify  our  notation  let  us  introduce  a  length  pQ  and  an  angular 
frequency  defined  by 


P 


o 


oj  =  8caa 
o  s 


[25] 


8 
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and  introduce  6t  =  a  6z.  In  terras  of  these  quantities,  we  can  easily 
s  s 

obtain 


5zP(p) 


Sz  a 

s 


Sts  exP|- 


expC-P2) 


[26] 


and 


a6z 


u)  2  -  2. 

—  P06ts  =  xpo6ts 
“o 


[27] 


In  eqs.  26  and  27  we  have  set  p”  *  p/pQ  and  X  =  w/«Q.  When  the 
scattering  phase  function  can  be  represented  by  eq.  23,  the  mutual 
coherence  function  y(w,p",ts)  will  satisfy  the  equation 


y(u,p,ts+6ts)  = 


2x«t, 


00 

f 

6ts  ) 

1  dp’p'  exp 
0 

Ti 

-2, 


r.  — 2  — ,2' 

exp  -  r1  “  *p 

L4  X«TS 

|J»fas)T(“-P,Ts) 

00 

=  2iy  J"  dp  'p  •  exp 
o 

— -  lexpC-p-2)  +  exp(-p"'2)| 

r  _2  —  2 

exp  -iy  (p  +  p’ 

)j  J0(2ypT’)  yGopP-,  ts) 

[28] 
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with  u  =  (4x<Sts) 

Choosing  the  step  size  6Tg  to  be  sufficiently  small,  one  can  use  the 
approximation 


exp 


fi 


exp(-p2)  +  exp( 


-P,2)| 


[• 


6ts  I  -2  -  2 

=  1  +  -  exp(-p  )+exp(-p'  )|[29] 


* 


in  eq.  28  such  that 


y(o),p,ts+6ts) 


=  2iu  /  dp'p'  jl  +  — ^ 


—2  —  2 
exp(-p  )+exp(-p  '  ) 


]! 


[30] 


exp^-iyfp2  p'2)^  JQ(2pp  p')  y(u),p,ts) 


Equation  30  is  the  relation  which  enables  us  to  generate  y(“»P» 
from  y(<*),F,ts). 

At  this  point,  it  may  be  worthwhile  to  comment  on  the  approxi¬ 
mation  in  eq.  29  in  comparison  to  the  quadratic  approximation 
2  2 

P(p)  ^  a  [1  -  k  p  /4a]  used  in  Refs.  5,  17.  In  eq.  29  it  is  evident 
5  #  2 

that  the  error  involved  is  of  the  order  5t  and  this  can  be  controlled 

s 

by  a  judicious  choice  of  the  step  size.  The  basis  of  using  the  quadratic 
approximation,  on  the  other  hand,  has  not  been  clarified.  In  fact,  the 
use  of  the  quadratic  approximation  in  the  present  treatment  will  lead 
to  incorrect  results. 
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4.0  STEPWISE  GENERATION  OF  THE  MUTUAL  COHERENCE  FUNCTION 

Based  on  eq.  30  derived  in  the  previous  section,  a  simple  alge¬ 
braic  relation  can  be  obtained  for  the  generation  of  the  mutual  coher¬ 
ence  function.  Symbolically  we  can  express  eq.  30  as 


y(w,p,ts+6ts)  =  [Gj  +  G2]  y(w.P.ts) 


[31] 


in  which  the  integral  transformation  operators  G^,  G2  are  defined  as  the 
following: 


00 

G1[y(u>,p’,ts)]  ■  2iy  J  dfT'fT'  exp[-iy  (p"2+p  ’2) ]  JQ(2yp  p  ')y(w,P  '  ,ts) 


[32] 


and 


St  r _  _2  _ 2 

G2[y(«,p’,ts)]  =  2iy  —  J  dp  V  [exp(-p  )  +  exp(-p '  )] 


[33] 


•  expl-iyfp^+p ,2)]  J  (2yp  p')  y(u,p',t  ) 


We  note  that  Gxy  and  G2y  are  of  the  zeroth  and  first  orders  of  6ts 
respectively. 

Consider  the  effects  of  applying  the  operators  Gj^  and  G2  to  the 
function  C  exp  (-ap2).  It  can  be  seen  that 
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GjtCexpf-oT2)] 


exp[-iy (p2+p ,2)]JQ(2yp  p ')Cexp(-ap  »2) 

[34] 


iyC 

(a+iy) 


exp 


Equation  34  can  also  be  expressed  as 


GjfCexpf-ap2)]  =  C'exp(-a'p2) 


[35] 


where 


iyC 

a+iy 


a  +  iy 


iy 


[36] 


Separating  C'  and  a'  into  real  and  imaginary  parts  such  that 
C'  -  CJ.  *  iC!  and  a'  ■  a^  +  ia!  ,  one  gets 


C'  = 
r 


yKaj+yJCy  -  arC.] 

2  '  7T~ 

\  ♦  (ou  +  y) 


C * 
Li 


u[arCr  +  (ai+y)Ci] 


T  , - ,2 

<*r  +  (Clj^+U) 


[37a] 


[37b] 
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2 

V  <* 


a'  = 


r  2  ,  .2 

°r  +  (Oj+tO 


“r  =  w  *  T 


M  (oij+y) 


“r  +  ^ai+w^ 


Similarly,  we  can  write 


G,[Cexp(-ap2)]  s  C ' 'exp (-a’ '^) 


[37c] 


[37d] 


[38] 


where 


C' ' 
r 


r  >  i 
l 


u[(g.*u)Cr  -  (l*»r)C1] 
(l+ar)2  +  (c^+v)2 

F[(l*«r)Cr  +  (oi+v)Ci] 
(l+ar)2  +  (cu+y)2 


V2(l+aJ 

"  =  1  ♦  - - j 

T  Cl+or)  ♦(Oj+y) 

2 

y  (a.+v) 

i  ^  2  2 

(l+ar)%(ai+y)^ 


[39a] 


[39b] 


[39c] 


[39d] 
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For  an  optical  path  of  thickness  Tg  =  n6tg,  the  mutual  coherence 
function  y(«,p\ts)  can  be  generated  from  the  relation 


Y(w,P,xs)  =  [Gj  +  G2]n  yG*>,P,0) 


[40] 


with  the  boundary  condition  y(u»P»0)  =1.  As  the  operators  Gj,  do 

not  commute,  the  expansion  of  [Gj  +  G2]n  gives  rise  to  2n  monomials  in 

G.,  and  G».  From  eqs.  35  and  38,  each  of  these  monomials  will  contribute 

a  term  of  the  form  C  exp(-ap  )  to  y(“»P>ts)  and  all  the  intermediate 

results  are  of  the  same  functional  form.  The  evaluation  of  these  terms 

can  now  be  carried  out  by  means  of  the  simple  algebraic  relations  in 

eqs.  35-39.  It  may  also  be  noted  that  in  the  expansion  of  [G^  +  G2]n 

there  are  nC  terms  of  the  m-th  order  in  6t  .  The  point  at  which  the 
m  s 

finite  series  can  be  terminated  in  practical  computation  is  effectively 
controlled  by  the  quantity  nCm(6Tg)m. 

For  large  values  of  n  the  quantity  nCm  increases  rapidly  with  m 

before  it  attains  its  maximum  value.  For  example,  if  n  =  30,  m  =  3  we 

30 

have  4060  terms  in  the  expansion  of  [G,  +  G~]  which  are  of  the  order 

3  A  ^ 

(6xg)  .  We  have  found  that  an  efficient  and  accurate  computational 

device  can  be  used.  This  consists  of  randomly  selecting  a  subset  of 

the  set  of  nC  terms  of  the  same  order  in  6x  for  m>3,  and  calculating 
m  s~ 

the  mean  value  of  their  contributions  to  yO»,p  =  0.  xg) .  The  total 

contribution  of  all  the  terms  of  order  (<5t  )m  is  then  approximated  by 
n 

C  times  the  mean  value  of  the  random  set. 
m 
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5.0  PULSE  PROPAGATION  IN  RADIATIONAL  FOG 

In  this  section  we  apply  the  previous  results  to  the  propagation 
of  a  laser  pulse  with  a  wavelength  equal  to  1.06  ym  in  radiational  fog. 

3 

It  has  been  shown  (Refs.  19,  20)  that  for  a  radiational  fog  with  0.1  g/m 
liquid  water  content,  the  scattering  phase  function  can  be  represented 
by  a  Gaussian  function 

v  2  -2  -1 

p(0)  =  exp(-ae^)  with  a  =  795.2  rad  ,  a  =  21.775  km  1 

and  6  denoting  the  scattering  angle.  Neglecting  the  small  absorptance. 
Figs. 1-3  show  the  results  for  the  mutual  coherence  functions 
y(w,F  =  0,Tg)  with  rg  =  2,  4  and  6  respectively  as  a  function  of  the 
frequency  x  3  u/uQ-  The  step  size  used  in  6Tg  was  0.2. 

Assuming  an  input  pulse  to  be  of  the  form 


Ii(t) 


exp(-t2/2T2) 


[41] 


it  can  be  readily  seen  that  the  output  pulse  is  given  by 


I 


IT 


00 


dto  exp 


cos  o>t  +  Yj  sin  cut] 


[42] 


where  yd  “d  YT  are  rea^  imaginary  parts  of  y(“>,p  =  0,  t  *  a  z) . 
Figures  4-6  show  the  output  pulses  for  the  cases  where  T  =  20  ns,  2  ns  and 
0.2  ns  with  all  the  pulses  normalized  in  such  a  way  that  max  I(t)  =  1. 

As  one  would  expect,  an  increase  in  the  optical  depth  results  in  an 
increase  in  the  pulse  distortion.  For  a  given  optical  depth,  on  the 
other  hand,  the  shape  distortion  in  shorter  pulses  is  more  severe. 
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TM  (to  Mita  m  r> 


FIGURE  5  -  Output  pulse  intensity  for  t 


TIME  («  wrili  •»  T  ) 


FIGURE  6  -  Output  pulse  intensity  for  t 
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6.0  CONCLUSION 


We  have  discussed  a  new  method  of  solving  the  parabolic  differen¬ 
tial  equation  of  the  mutual  coherence  function  based  on  the  split-step 
Fourier  method.  The  principal  result  is  a  set  of  simple  algebraic 
recurrence  relations  from  which  the  mutual  coherence  function  can  be 
generated  for  an  arbitrary  optical  depth.  Compared  to  previous  work, 
the  present  treatment  does  not  require  any  additional  assumption  to 
solve  the  parabolic  equation.  The  results  are  applicable  to  turbid 
media  provided  the  typical  scatterer  dimension  is  comparable  or  larger 
than  the  pulse  wavelength. 
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numerical  algorithms  in  the  form  of  algebraic  recurrence  relations  are  numerical  algorithms  in  the  form  of  algebraic  recurrence  relations  are 
derived  for  the  study  of  the  phenomenon  of  pulse  broadening  induced  by  derived  for  the  study  of  the  phenomenon  of  pulse  broadening  induced  by 
multiple  scattering.  (II)  multiple  scattering.  (II) 


